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Summary 

The  motion  of  an  edge  crack  extending  non-uniformly  in  an  elastic  half¬ 
space  under  conditions  of  anti-plane  shear  is  analyzed.  An  expression  for 
the  stress  Intensity  factor  at  the  crack  tip  is  obtained,  and  an  energy 
balance  crack  propagation  criterion  is  used  to  find  the  eguation  of  motion 
of  the  tip.  On  solving  this  equation  numerically,  it  Is  found  that  crack 
arrest  occurs  before  the  second  reflected  wave  from  the  boundary  reaches 
the  tip. 

In  the  second  half  of  this  Investigation,  a  numerical  procedure  for 
studying  anti -plane  shear  crack  propagation  problems  using  finite  differences 
Is  developed.  To  approximate  the  elastodynamic  field  as  accurately  as  pos¬ 
sible  near  the  moving  crack  tip,  where  singular  stresses  occur,  the  local 
asymptotic  displacement  field  near  the  tip  is  Incorporated  into  the  finite 
difference  scheme.  The  numerical  procedure  is  applied  to  the  edae  crack 
problem  analyzed  In  the  first  part  of  this  study,  and  the  numerical  and 
exact  results  are  compared. 


The  results  communicated  In  this  paper  were  obtained  in  the  course  of  an 
Investigation  supported  In  part  by  Contract  N00014-75-C-0196  with  the  Office 
of  Naval  Research  in  Washington,  O.C. 
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GENERAL  INTRODUCTION 


The  analysis  of  rapid  crack  propagation  and  crack  arrest  in  solid 
bodies  has  recently  been  generating  a  considerable  amount  of  interest. 
Generally  the  failure  of  a  structure  through  fracture  is  avoided  by  pre¬ 
venting  the  onset  of  unstable  crack  growth.  This  can  be  done  by  ensuring 
that  the  stresses  do  not  exceed  certain  material -dependent  limits  imposed 
by  critical  flaw  sizes  in  the  material.  However  excluding  unstable  crack 
extension  under  a  wide  range  of  operating  conditions  can  often  be  very 
difficult  and  costly.  A  suitable  second  measure  of  safety  would  be  to 
ensure  the  timely  arrest  of  a  propagating  crack  and  this  is  often  a 
practical  necessity  in  important  structures  such  as  nuclear  pressure 
vessels,  welded  ships  and  pressurized  pipelines.  In  cases  like  these, 
unstable  crack  propagation  without  arrest  could  have  very  serious  con¬ 
sequences  . 

Dynamic  crack  propagation  problems  generally  fall  into  two  categories. 
Firstly,  for  prescribed  crack  history,  the  elastodynamic  field  in  the  body 
is  to  be  determined.  Alternatively,  for  an  unstable  crack  in  a  body, 
given  all  the  necessary  Information  concerning  the  geometry,  the  loads 
and  the  material  of  which  the  body  is  composed,  together  with  a  crack 
propagation  criterion  as  well,  one  can  attempt  to  determine  the  subsequent 
motion  of  the  crack.  The  latter  case  will  concern  us  in  this  study. 

Now  the  rapid  motion  of  a  crack  in  a  body  introduces  inertia  effects 
in  the  form  of  stress  waves  radiating  outward  from  the  moving  crack  tip. 
Hence  in  analyzing  a  dynamic  fracture  problem,  the  reflections  and  dif¬ 
fractions  of  these  stress  waves  by  the  material  boundaries  and  by  the 
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crack  Itself  must  be  taken  into  account.  Because  of  the  inherent  mathe¬ 
matical  difficulty  in  the  boundary-initial  value  problem  due  to  the  presence 
of  the  propagating  crack  and  the  disturbances  caused  by  these  wave  motions, 
the  theoretical  analysis  of  such  problems  is  exceedingly  difficult  and  no 
generally  applicable  solutions  are  available,  except  for  some  simple 
regions.  Section  1  Includes  a  discussion  of  some  of  the  analytical  work 
which  has  been  done. 

Consequently,  in  view  of  the  considerable  interest  in  dynamic  frac¬ 
ture  problems,  numerical  treatments  of  these  problems  are  of  prime  impor¬ 
tance.  It  Is  essential  that  reliable  numerical  procedures  be  developed 
which  are  relatively  easy  to  use  and  are  capable  of  accurately  simulating 
dynamic  crack  propagation.  This  is  by  no  means  an  easy  task  in  view  of 
the  moving  crack  tip  and  the  singularity  in  the  stress  and  velocity  fields 
at  the  tip  Itself.  Section  6  considers  some  of  the  numerical  techniques 
used  In  studying  crack  propagation  problems. 

The  work  contained  In  this  Investigation  is  restricted  to  infini¬ 
tesimal  deformations  of  homogeneous  Isotropic  elastic  bodies  and  to  anti¬ 
plane  shear  displacement  fields.  Although  of  less  practical  Interest, 
anti-plane  shear  Is  fundamentally  simpler  both  analytically  and  numerically 
than  plane  strain,  having  just  one  non-vanishing  displacement  component 
and  one  characteristic  wave  speed.  Since  certain  basic  features  of  both 
plane  and  anti-plane  crack  propagation  problems  are  similar,  the  anti-plane 
shear  problems  do  yield  considerable  insight  Into  crack  propagation  in 
general  and  may  be  used  as  a  qualitative  guide  In  analyzing  the  more  dif¬ 
ficult  plane  strain  problems. 

In  Part  I  of  this  study,  the  motion  of  a  propagating  edge  crack  in 
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an  elastic  half-space  under  conditions  of  anti-plane  shear  is  analyzed. 

An  energy- balance  propagation  criterion  is  used  at  the  crack  tip  to  de¬ 
termine  the  velocity  of  the  tip  and  an  expression  for  the  stress  intensity 
factor  at  the  moving  crack  tip  is  obtained. 

A  numerical  procedure  for  analyzing  anti -plane  shear  crack  propaga¬ 
tion  problems  using  finite  differences  is  developed  in  Part  II.  The  local 
asymptotic  displacement  field  near  the  moving  crack  tip  is  incorporated 
into  the  finite  difference  scheme  so  that  in  the  vicinity  of  the  crack 
tip,  where  singular  stresses  occur,  the  elastodynamic  field  is  approximated 
in  the  numerical  process  as  accurately  as  possible.  The  numerical  scheme 
is  applied  to  the  edge  crack  problem  studied  in  Part  I  and  the  numerical 
and  exact  results  are  compared. 


I.  THE  ANALYTICAL  SOLUTION  OF  THE  EDGE  CRACK  PROBLEM 


1 .  Introduction 

The  existing  theoretical  solutions  to  dynamic  crack  propagation  problems 
have  limited  applicability,  since  they  rarely  include  the  effects  of 
external  boundaries.  For  non-uniform  crack  velocity,  exact  solutions  of 
dynamic  crack  propagation  problems -val id  for  all  time-have  only  been 
found  for  simple  configurations  such  as  that  of  a  semi -infinite  crack  in 
an  otherwise  unbounded  body.  Generally,  constant  velocity  crack  problems 
are  more  amenable  to  analysis  and  have  been  studied  using  a  number  of 
techniques.  However,  to  determine  mathematically  how  a  plane  crack  propa¬ 
gates  in  a  given  situation,  first  the  elastodynamic  field  near  the 
crack  tip  has  to  be  found  for  the  tip  moving  in  an  arbitrary  fashion  in 
its  plane  and  then  a  crack  propagation  criterion  must  be  employed  to 
find  the  actual  motion  of  the  tip. 

The  linear  elastodynamic  field  of  a  propagating  crack,  extending  non¬ 
uni  formly  under  general  anti -plane  loading,  was  investigated  by  Kostrov 
[1]  and  Eshelby  [2],  while  Freund  [3]  and  Kostrov  [4]  analyzed  the  plane 
strain  case.  Kostrov's  procedure  in  the  anti-plane  problem  [1]  has  since 
been  used  successfully  by  a  number  of  authors  and  Is  employed  here. 

Kennedy  and  Achenbach  [5]  consider  the  extension  of  a  crack  which  was 
formed  by  cutting  at  a  high  velocity  into  the  surface  of  an  elastic  solid 
subjected  to  a  uniform  anti-plane  loading.  Dynamic  crack  growth  generated 
by  an  anti-plane  shear  wave  incident  on  a  crack  tip  was  investigated  by 
Achenbach  [6,7].  Molchanov  and  Nikitin  [8],  considering  a  finite  crack 
propagating  in  an  unbounded  body  subjected  to  anti-plane  shear,  studied 
the  effect  of  the  stress  waves  emitted  from  one  of  the  moving  crack  tips 


upon  the  motion  of  the  other  tip,  just  after  the  arrival  of  the  first 
such  wave.  Recently  Melville  [9]  examined  approximately  the  effect  of 
reflected  stress  waves  on  the  stress  intensity  factor  for  an  edge  crack 
in  anti-plane  shear  during  crack  propagation  and  after  crack  arrest.  In 
simplifying  the  analysis,  he  assumed  a  constant  crack  tip  velocity  as 
well  as  a  simple  form  for  the  stress  distribution  ahead  of  the  initial 
crack  tip. 
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Consider  an  elastic  half-space  containing  a  traction- free  edge  crack 
of  initial  length  4  along  the  x-axis  (Fig.  1(a)).  The  deformation  is 
assumed  to  be  one  of  anti-plane  shear  and  hence  the  only  non-zero  displace¬ 
ment  component,  which  is  denoted  by  w,  is  in  the  z-direction  out  of  the 
plane  of  the  paper  in  Fig.  1(a).  Initially  the  body  is  in  a  deformed 
equilibrium  state,  with  the  crack  tip  held  fixed  at  x  =  2,Q.  The  initial 
Static  displacement  field  of  the  body  in  this  state  satisfies  the  boundary 

condition  on  x*0  given  by  w  =  w„  for  y  >  0  and  w  * -w  for  y  <  0, 

o  o 

where  wQ  is  assumed  constant.  This  fixed  boundary  condition  at  x*0 
Is  maintained  for  all  time.  At  time  ts0  the  crack  tip  is  released  and 
crack  propagates  along  the  x-axis,  assuming  wQ  is  sufficiently  large. 

The  position  of  the  crack  tip  is  given  by  xa4(t)  at  time  t. 

Since  w«w(x,y,t),  the  stress  tensor  has  only  two  non-zero  components 
given  In  terms  of  the  displacement  gradients  by 

°ZW  *  UW,X  3nd  azy  *  yW,y  ’  (2-1) 

y  being  the  Infinitesimal  shear  modulus  of  the  material.  The  displacement 
equations  of  motion  for  a  homogeneous  Isotropic  linear  elastic  body  reduce 
to  the  wave  equation 

"•»x*w>yy  ’  cZ'v/B  ’  (2'2) 

c 

where  p  Is  the  mass  density  and  c  is  the  shear  wave  speed  in  the  body. 
The  crack  tip  position  4(t)  Is  not  specified  a  priori  but  is  later  de¬ 
termined  using  an  energy- balance  rropagatir  criterion.  The  only  restric¬ 
tion  Imposed  on  4(t)  Is  that  the  c^ck  tip  velocity  be  less  than  the 
shear  wave  speed  l.e.  0<i(t)  <c. 
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During  crack  propagation,  the  moving  crack  tip  acts  as  an  energy  sink. 
As  new  fracture  surface  is  being  formed,  the  loss  of  mechanical  energy  re¬ 
sults  from  the  release  of  the  internal  cohesive  bonds  between  adjoining 
particles  in  the  material.  The  energy  release  rate  G  is  defined  as  the 
amount  of  energy  being  absorbed  by  the  moving  crack  tip,  per  unit  crack 
tip  velocity  e(t)  and  per  unit  thickness  of  the  crack.  In  general  G 
depends  upon  t,  i(t)  and  the  crack  history  jj.( r ) ,  O^rnt.  The  flux 
of  energy  into  the  tip  per  unit  time  is 


F  =  G*  ,  i*t(  t) 


From  energy  balance  considerations,  one  has 


(2.3) 


F.P-£(W+K)  . 


(2.4) 


where  W  and  K  are  the  strain  and  kinetic  energies  of  the  whole  body 
respectively,  and  P  represents  the  rate  at  which  work  is  being  done  by  the 
external  forces  on  the  body.  Atkinson  and  Eshel by  [10] showed 

F'USf  flwlTa[7((tr)2*(l7)2)+2(l7,2>x}dS  •  (2’5) 


where  S  is  a  closed  contour  surrounding  the  crack  tip  and  moving  with  it, 

and  n=(n  ,  n  )  is  the  unit  outward  normal  to  S. 
y 

For  a  crack  along  the  x  axis  with  tip  at  x=  t(t)  ,  the  dynamic 

anti-plane  shear  stress  intensity  factor  Kd=Kd(t)  can  be  defined  through 


°Zy(x*0,t)  =  (2Tr)"\(tHx  -  +0(x  -  i(t))H 


(2.6) 


as  x-*t(t)  .  Writing  the  asymptotic  expressions  for  the  stress  and  d i s ■ 


I*  -]  -J 


placement  fields  near  the  crack  tip  in  terms  of  ,  equation  (2.5)  with 
(2.3)  gives  (see  [2J) 


(2.7) 


G  as  given  by  (2.7)  will  be  employed  in  formulating  the  crack  propagation 
criterion  to  be  used  for  the  determination  of  the  crack  history  i{ t)  . 
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3.  Solution  of  the  static  problem 

To  investigate  the  crack  tip  motion,  the  initial  static  stress 
field  ahead  of  the  crack  is  required.  The  configuration  and  boundary 
conditions  are  shown  in  Fig.  1(a).  The  crack  faces  are  traction  free, 
so  that 


w,y  =  0  on  y  =  0±  ,  0  ^x<iQ  .  (3.1) 

Since  the  displacement  w  is  now  time- independent,  (2.2)  reduces  to 
Laplace's  equation  which  must  be  solved  in  the  right  half-plane  excluding 
the  crack  and  subject  to  the  appropriate  boundary  conditions.  It  is  also 
required  that  w  be  bounded  near  the  crack  tip. 

The  mapping 

t*  arcsin(l  -  z2/Zq)$s  (3.2) 

z  =  x  +  iy  ,  t  =s  +  in  ,  (i  =  /-T) 

2  2  V- 

on  taking  the  branch  cut  for  (1-z  /iQ)a  along  the  positive  x  axis, 
transforms  the  configuration  in  Fig.  1(a)  onto  the  semi-infinite  strip 
in  Fig.  1(b)  with  the  boundary  conditions  as  shown.  It  can  now  be 
readily  verified  that  the  static  displacement  field  wst=  wst(x,y)  is 


.st 


2w. 


2m. 


w 


Ret 


Re 


arcsin( 1  -  z2/d2)^ 


(3.3) 
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The  stresses  are  found  by  differentiation  to  be  given  by 


st  ,  st  =!To  n._2/52r4 
°zx  1azy  Treo  z  'V 


(3.4) 


Ahead  of  the  crack  tip,  we  have 


I*(*.0)-»|^l(x,0)-^(x2-^rls.  x>i  .  (3.5) 


For  an  equilibrium  crack  along  the  x  axis  with  tip  at  x =  iQ  , 
the  static  stress  intensity  factor  K  in  anti-plane  shear  is  defined 


through 


o“(x,0)  =  (2.)  J»KJt(x-  l0)'H  *0(x  -  l0)H 


(3.6) 


as  x-*-zQ  ,  and  the  static  energy  release  rate  G$t  is  given  by  (see[ll]) 


Gst  *  Yt^ 


Then  using  (3.5)  for  the  edge  crack  under  consideration,  we  have 


K  = _ a_ 

st  U*0)'5 


(3.7) 


(3.8) 


St  irl. 


(3.9) 


•  •  *  •  -  *  *  %  •  . *■„  A. 
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A.  Formulation 


The  moving  crack  tip  acts  as  a  source  of  stress  waves  which  propagate 

2  2  2 

outward  from  the  tip  behind  the  cylindrical  wavefront  (x-£q)  +y  =(ct)  . 
Using  s  =  ct  in  (2.2),  where  c  is  the  shear  wave  speed,  the  dynamic 
displacement  field  w(x,y,s)  resulting  from  this  disturbance  satisfies 


Let  T(s)  denote  the  position  of  the  crack  tip  in  terms  of  s  at  time  t=s/c, 
so  that  T(s)  *  £(s/c ) .  In  T(s),  wenowdrop  the  bar  for  convenience  and  from 
here  on,  s  will  be  used  as  a  normalized  time  variable.  The  boundary  conditions  are 


w,y  =  0 


0  sx<n(s) 


(4.2) 


wQ  ;  x  =  0  ,  y  >0 
-w  ;  x  =  0  ,  y  <  0 


(4.3) 


Initially  at  s  =  0  ,  the  body  is  at  rest  i.e.  w,s(x,y,0)=0  and  the 
initial  static  displacement  w  =  wst(x,y)  is  given  in  (3.3). 

Let  w'(x,y,s)  »  w(x,y,s)  -  wst(x,y)  denote  the  change  in  displace¬ 
ment  from  the  initial  static  configuration.  The  boundary  condition  at 
x*0  now  becomes  w'(0,y,s)  =  0.  As  a  consequence  of  this  fixed  end  con¬ 
dition  at  x*0,  the  reflected  stress  waves  must  be  equal  in  magnitude  but 
opposite  in  sign  to  the  incident  waves  at  this  location.  To  incorporate 
the  effect  of  these  reflected  waves,  we  consider  a  center  crack  in  an 
infinite  body  with  tips  at  x*±l{ s)  i.e.  moving  at  the  same  rate  but  in 


opposite  directions.  Then  the  residual  displacement  w'  must  satisfy  the 
two-dimensional  wave  equation  (4.1)  In  the  full  (x,y)  -  plane  excluding  the 
crack,  with  zero  Initial  conditions  and  the  boundary  conditions  on  the 
crack  surfaces  are  given  by 

0 ,  0  a  x  <  l 

0  ,  (4.4) 

-P(x),  lQ  <  x  < £(s) 


u  |f  u.oM 


>■  % 


0,  -£q<X£Q 

P(-x),  -£(  s)<x<-£„ 


(4.5) 


where, from  equation  (3.5),  P(x)  is  given  by 


P(x)*o^(x,0)  = 


2yw 


o  /  2  „2 


-i 


*  <*' •  x>to 


(4.6) 


The  solution  satisfies  w'(0,y,s)  =  0  as  required.  Hence  the  half-plane 
problem  with  the  edge  crack  under  consideration  is  equivalent  to  the 
full -plane  problem  with  a  center  crack  as  described  above. 


B.  Green's  Function  method 

Since  the  initial  static  displacement  field  w(x,y)  is  anti -symmetric 
in  y,  w(x,y,s)  and  hence  w'(x,y,s)  is  also  anti -symmetric  in  y.  In 
what  follows,  only  the  half-plane  y^O  is  considered.  Kostrov's  method 
of  solution  [1],  employed  here,  uses  the  Green's  Function  G  which  satisfies 
the  two-dimensional  wave  equation  (4.1)  in  the  half-plane  y^O,  with  zero 
initial  conditions  and  the  boundary  condition 
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<£5i 

h\ 


m 


C5M 


ffzy(x,0,s)  =  y  —  (x,0,s)  =  6(x  -  xq)6(s  -  sQ)  ,  (4. 


where  -co<x0<oo  and  OssQ<s  .  The  displacement  wave  generated  by  this 


impulsive  load  is  (see  [13]) 


G(x-x0,y,s-s0)  -  -  H{(s  -  s0)  -  [(x  -  xQ)2  +  y2]^J  , 


where  H{  }  is  the  Heaviside  step  function  and  R  is  given  by 


This  is  a  cylindrical  wave  emanating  from  the  line  x  *  xQ  ,  y  =  0  for 


(4.8) 


(4.9) 


s>so  ‘ 


For  general  boundary  conditions  of  the  form 


°Zy(XtO,s) a ct(x,s)  , 


(4.10) 


linear  superposition  gives 


,  rr  a(x_, s) 

w' (x,y,$)  *  -  —  ||  - r2_dxodso  ’ 


(4.11) 


where  R  is  given  in  (4.9)  and  the  region  of  integration  A  is  the 
parabolic  region  in  the  xQ  -  sQ  plane  satisfying 


s-so-[(x-xo)2+y2]!j20  ,  s>sQiO  . 


(4.12) 


On  y*0  ,  A  reduces  to  a  triangular  region.  Introducing  the  character¬ 
istic  coordinates 


5=(s-x)/»/5’  and  n=(s  +  x)/>/2"  ,  (4.13) 

the  expression  for  R  in  (4.9),  on  setting  ya0,  simplifies  to 

R*^(5-C0))5(n-n0)i5  .  (4.14) 


From  (4.11),  the  displacement  w'(x,0,s)  a w' ( C ,n)  ,  on  using  (4. 14),  becomes 


w'U  n)  * 


1 

iry/2" 


•  <*Q  f  g^0»^^Q 

<e-«0)H  >  ln-nort 


(4.15] 


The  (x.s)  -diagram  is  shown  in  Fig.  2,  along  with  the  crack  tip  trajectories 
and  wavefronts.  To  determine  the  motion  of  the  crack  tip,  the  energy  re¬ 
lease  rate  and  hence  the  stress  Intensity  factor  at  the  tip  must  be  found. 

C.  The  field  before  the  arrival  of  the  first  reflected  wave 

The  first  reflected  stress  wave  strikes  the  crack  tip  at  s » sp  , 

where  s„  satisfies 
r 

*(sr)  +  *0*s  r  .  (4.16) 


From  [2],  or  following  the  approach  used  in  [12],  the  shear  stress  ahead 
of  the  crack  tip  before  the  arrival  of  the  first  reflected  stress  wave  is 


where 


;•  % 


Us) 


■o‘S> 


P(v)dv 


(Us)  -  v)5 


(4.18) 


Evaluating  IQ  ,  where  P(v)  *  cr^(v,0)  is  obtained  from  (3.5),  gives 


4uw 


°  mi  k 

L  >  ^ 

U+  iQ) 2 


i  -  i. 


H 


i  +  i. 


(4.19) 


Here  K(k)  is  the  complete  elliptic  integral  of  the  first  kind.  The  stress 
intensity  factor  then  is 


Kd(s)»(2/1r)ls(l-U(s))i5I0(s)  ,  Oss<S| 


(4.20) 


and  the  energy  release  rate  at  the  crack  tip  is 


•  «»<», 


(4.21) 


D.  The  field  after  the  arrival  of  the  first  reflected  wave 

After  the  first  reflected  stress  wave  strikes  the  moving  crack  tip  at 


the  Instant  s  ■  sr  given  by  (4.16),  a  diffracted  wave  is  emitted.  This 


diffracted  wave  propagates  outward  behind  a  cylindrical  wavefront  with 
the  shear  wave  speed,  hits  the  fixed  boundary  at  x  =  0  and  is  reflected. 

It  strikes  the  crack  tip  again  at  s s sd  ,  where  from  Fig.  2,  sd  satisfies 


sd't(sd)"sr  +  t(sr)  ' 


(4.22) 


! 


.*«•  *J<*  . 

■  .  r,  •  .  ♦  . 
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We  now  analyze  the  stress  field  in  the  crack  plane  near  the  tip  for 


sr<s<sd  . 


In  terms  of  5  and  n  .  the  crack  tip  trajectory  x=n(s)  becomes 
n  *  L( c )  where  L( c)  satisfies 


UShiS.  =  z/Uil 


1  {  =**'  +  K\ 

n  \  rr  ) 


(4.23) 


Similarly  x*-fc(s)  becomes  n  =  LU)  where 


Ui. LlL  a  -t  l  Lj-O-ti  ^ 


ft 


\  G 


'4.24) 


The  initial  crack  tip  locations  x  =  ±£  are  denoted  by  nsn+(c)s5^  £Q 
Consider  the  point  with  characteristic  coordinates  (s2»  n2)  shown 


In  FTg.  2,  where  for  ?r*t0/i^T  , 


cr<52<5d  and  n2>L(52) 


(4.25) 


Here  Is  the  value  of  5  when  the  diffracted  wave  meets 

Q  U  Q 


the  right  crack  tip  at  s  *  sd  and  x=£(s{j).  The  domain  of  dependence  is  be¬ 


low  the  dotted  1  ines  shown.  Wedenoteby  a^C.n)  the  shear  stress  ozy  at 


the  point  U,n)  ahead  of  the  right  crack  tip  for  -C  <£<C  and  similarly 


let  o^(c,n)  represent  the  shear  stress  aZy  ahead  of  the  left  crack  tip 


for  -nr<n<nr  ;  where  nf  *  i0/»^  .  The  change  in  displacement  w'U2,n2) 
ahead  of  the  right  crack  tip  is  zero  giving,  from  (4.15) 


t  >y  v-v.'  •  :<  .  v •••  -wo  %• •  vivVvv-:  • ;  ;  . 


dnd?  -  0 


(4.26) 


C2  n_(s) 

P(S,n) 

5r  L(?) 


Us) 


R 


Here  R  is  given  in  (4.14), 

P(S,n)*P((n-S)/^)  =  P(x) 

and  P(s.n)  -  P((n  -  s)/#)  ■  P(x)  »  -P(-x) 


(4.27) 

(4.28) 


At  (S-|  ,n-j )  in  Fig.  2,  where  -Sr<S1<Sr  and  n-,>L(s1)  ,  that 
is  before  the  arrival  of  the  first  reflected  wave,  w'(Si»n-|)  is  also 
zero,  giving 


*1  L<«>  51  ”1  > 

f  {  ^^-dnds-f  |  °]  ft’  dnds  *  0  .  (4.29) 

-Sr  n+(s)  -Sr  Us) 

Extending  the  definition  of  a-|(s»n)  to  the  region  sr<S<S2  and  n>L(s). 
where  o^(s.n)  is  here  taken  as  the  shear  stress  at  (s,n)  if  the  first 
reflected  stress  wave  is  not  taken  into  account,  then  the  second  term  of 
(4.26)  can  be  rewritten  as 


By  virtue  of  (4.29)  with  ( C i » n-j )  replaced  by  ( ^2 » ^2 )  snd  (4.30),  equation 
(4.26)  on  using  (4.14)  becomes 


£  n2 

r  f  Ca2(e»n)-ai(c»n)Ddn 

€r^2“5^  L(e)  (n2-n)2 


*2  £,  } 
f  — d6  |  P(C.n)dn  _  f  d£  f  gl(*>»^dri 


'  (52‘5) '  r  .  ( Vn)  £  (C2-«r  J  (n2- 
■r  £  L(£)  *  cr  2  -fY  2 


(4.31 


Equation  (4.31)  is  clearly  satisfied  by 


n2  n  (5) 

r  [<J2(c»n)-a1(£,n)]dn  7 

(n2-n)1*  J  (n? 

LU)  2  L(£)  2 


L(C) 

P(6.n)dn  .  f  alU>n)dn 


-  - - r—  •  (4.32) 

(n2-n)  J  (n2-n)  3 
r 


This  is  an  integral  equation  for  o2U,n)  -o^C.n)  of  the  Abel  kind  which 
can  be  solved  analytically  to  give 


aotS’n)  -OiU.n)  =  - - r 

2  1 


ruiunuu-u 
J  n  -  u 


S^.uJWO-u]1 


Similarly  from  (4.29)  with  (c-j.n-j)  replaced  by  Ug.T^) ,  it  follows  that 


f  0^5, n)  L|5)  p(  v 

— - irdn  =  — l  dn  , 

(n2‘n)a  (n?-n)2 

L(0  2  n+(s)  2 


which  upon  solving  gives 


0i(c»n)  =  t 

ir(n-L(c))is 


n+U) 


Combining  (4.33)  and  (4.35)  gives  the  change  in  shear  stress  ahead  of  the 
right  crack  tip,  taking  the  first  reflected  stress  wave  into  account  i.e. 


«,U.n) - ! — r  J  f  l.  du 

.(n-L(c)^  J(l) 


n  (5) 


diiru^-u]1* ...  _V°  51(e.“)[L(«)-«]‘‘ 


«"fc  W»"“  •  •  ******  •  •  »*•*•*■  .*'  iVAV' 
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s(x,s)  -  s  -  (1  -  a'(s))_1(x-Ji(s))  +0(x  -  *(s))2  ,  (4.42) 

whence  i(s)  *  £( s )  +o(x  -  £( s ) )  .  From  (4.41)  6-|(x,s)  =  F(s) +o(x  -  t(s) ) , 

where 

F(s)  -  -  (£(s)-s-a0)/2  .  (4.43) 

Equation  (4.39)  may  be  rewritten  as 

S2  + t(?2)  *  s  -  £(s)  -  (x  -  t(s))  .  (4.44) 

Let  s2  be  the  solution  of 

s2  +  t(s2)  s s  -  t(s)  .  (4.45) 

For  srss<s£j,  one  can  show  that  Oss2<sr.  From (4.44) ,  ?2=s2+0(x-£(s) )  and 
so  i(s g)  *  t(s2) +0(x-t(s)) .  Physically,  the  stress  wave  arriving  at  the 
right  crack  tip  xsi(s)  at  time  t=s/c  was  emitted  at  the  left  crack 
tip  x=-*(s2)  at  time  t2  =  s2/c  .  Using  the  above  asymptotic  results 
for  s’  and  ?2  near  xsi(s),  equation  (4.40),  after  some  further 
manipulation,  becomes 

o,(x,s)  -  l(s)+d*-l(s))H  ,  (4.46) 

2  w(x-»(s))’!  r 


.-S  '.  >v I 
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where 


P(v) dv 
(t(s)-v)‘s 


P(v)dv 

(l(s)-v)h 


-£(s2) 


-6(s) 


( v,s-£(s)+v)dv 
U(s)-v)1* 


(4.47) 


Here  fF(s)  is  given  by  (4.43)  and  s2  satisfies  (4.45). 

In  (4.47),  a^(x,s)  the  shear  stress  ahead  of  the  left  crack  tip 

for  Oss<sr  is  required.  Equation  (4.35)  gives  c^U.n),  the  shear  stress 
ahead  of  the  right  crack  tip  before  the  first  reflected  wave  arrives. 
Transforming  a^s.n)  to  the  variables  x  and  s  as  was  done  for  c^U.n) 
gives 


<J-|  (x,s)  = 


1 


£(?) 


ir(x-i(s) ) 


PLy)[,^s)-y,lh  dv 

x  -  v  av  » 


(4.48) 


where  F  is  given  by  (4.38).  By  virtue  of  the  symmetry  in  the  problem, 

replacing  x  by  -x  ,  v  by  -v^  and  P(v)  by  -P(-v.j)  in  the 

integrand  of  (4.48)  gives  d^(x,s),  the  shear  stress  ahead  of  the  left 

crack  tip  for  Oss<s„  ,  i.e. 

r 


a^x.s) 


1 

w(-x-a(s1))>5 


). 

P(v1)[t(s1)+v]] 


>5 


(4.49) 


in  which  P( v-j )  s  -P(-v^)  and  s^  satisfies 


s1  -  t( s ^ )  =  s  +  x 


(4.50) 
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Equation  (4.49)  now  yields 


d1(v,s-ii(s)+v)  =  - 


tt(-V-£(s1  )  )' 


P(v1)[£(s1)+v1]i5 


■4(5,) 


dv,  ,  (4.51) 


where  from  (4.50),  is  obtained  from  the  relation 


s-,  -  fi,(s1 )  =  s  -  £(s)  +  2v 


(4.52) 


For  srss<sd  ,  one  can  show  Oss2<sr  and  Oss-|<s2  .  Evaluating  the 
integral  in  (4.51),  one  gets 


3,(v,s-£(s)+v)= - - J - - 

»  (£(s1)+£0)^(-v-£(s1))15 


(4.53) 


where 


.0.-0  *(§,)-£„  £(s,)-£ 

k  *«*i>  ‘Hl^  a"d 


(4.54) 


Here  n(a,lc)  is  the  complete  elliptic  integral  of  the  third  kind.  Note 
that  asO  for  all  v  in  the  interval  (-ff(s) ,  -  £(s2))  .  The  second  term 
of  (4.47),  on  performing  the  integration,  is 


-MS,) 


U(s)-v) 


17  U(s)n0p 


(4.55) 


In  (4.55)  F(<p,k)  is  the  incomplete  elliptic  integral  of  the  first  kind. 


• 

C  vV  ■ 


•\.s  »'*  •  .V.  . 
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5.  Motion  of  the  crack  tip 

During  the  fracture  process,  the  moving  crack  tip  absorbs  energy 
from  the  surrounding  material.  According  to  prevailing  views  of  fracture, 
the  energy  absorbed  by  the  tip  is  used  to  create  the  new  free  surface 
generated  by  crack  propagation.  In  the  usual  quantitative  model  of  this 
process,  one  postulates  the  existence  of  a  "specific  fracture  energy"  r 
which  is  characteristic  of  the  material  at  hand  and  which  has  the  property 
that,  during  crack  propagation,  the  relation 

G  =  r  (5.1) 


holds  at  every  instant.  According  to  (4.21),  (4.19)  and  (4.59),  (4.57), 

G  depends  on  the  crack  length  history  i(r) ,  Osrss  and  on  the  in¬ 
stantaneous  crack  tip  velocity  *.'($).  In  general  r  may  be  a  (given) 
function  of  «/(s);  we  consider  only  the  simplest  model,  in  which  r 
is  an  absolute  constant.  Equation  (5.1)  with  r s constant  is  thus  our  crack 
propagation  criterion.  It  may  be  regarded  as  an  "equation  of  motion"  for 
the  crack  tip. 

In  Section  2,  it  was  assumed  that  the  crack  tip  was  held  fixed  until 
time  ta0,  when  it  is  released  and  the  crack  propagation  process  begins. 
For  propagation  to  occur  at  all,  (5.1)  requires  that  the  static  energy 
release  rate  G$t  satisfy  G^r.  Thus  by  (3.9),  one  must  have 
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3 


a 


A  dimensionless  parameter  Q  is  now  introduced  through 


2yw‘ 

Q  =  — - 

^  —  n  v< 


(5.3) 


Then  for  crack  extension  to  occur,  from  (5.2)  it  is  necessary  that  Q 
be  larger  than  unity. 

Before  the  arrival  of  the  first  reflected  stress  wave,  the  crack  tip 
equation  of  motion  is,  by  (4.21),  (4.19)  and  (5.1), 


,6v 


l-!'(s)P  K(k] 


ln'(s )/  (T(I)n  ) 


(5.4) 


Using  (5.3)  and  the  definition  of  k*k(i(s))  in  (4.19),  the  differential 
equation  in  (5.4),  on  solving  for  i'(s),  becomes 


*'(s) 


,  Qss< 

K(kr(1-kz)z-hr4/16Q2  r 


(5.5) 


At  s*0,  we  have  k*0  and  K(0)  *  tt/2  giving 


t'(o)=9^1  . 

(T  +  l 


(5.6) 


Since  Q>1  is  necessary  for  crack  extension  to  occur,  t'(0)  must  be 
positive.  Because  k*kU(s)),  equation  (5.5)  is  an  autonomous  differential 
equation  for  t  with  a  critical  point  at  k = kr,  where 


Klk/fl-k2)  -,2/4Q  . 


(5.7) 
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Since  an  infinite  amount  of  time  is  required  to  reach  a  critical  point,  the 

crack  cannot  arrest  before  the  arrival  of  the  first  reflected  stress  wave. 

At  sssr  ,  for  Q>1  the  crack  tip  is  moving  with  non-zero  velocity 

*'(sj  where  0<t'(s  )<1  . 
r  r 

After  the  arrival  of  the  first  reflected  wave  at  s  =  s  until  the 

r 

arrival  of  the  second  reflected  wave  at  s = sd  ,  from  (4.59)  and  (5.1), 
the  equation  governing  the  crack  tip  motion  is 


Here  I(s)  is  given  by  (4.57).  Using  (5.3),  this  equation  reduces  to 
solving 


«'<S)  - 

I(s)Wl6QZ 


*s<sd 


(5.9) 


where  fc(sr)  is  obtained  from  the  solution  of  (5.5).  From  (4.57),  I(s) 
depends  upon  the  instantaneous  variables  s  and  i( s),  and  upon  the 
delayed  crack  length  ^(sg).  If  for  some  value  of  s  in  the  interval 


I(s)2  -  tt2/4Q  =  0 


(5.10) 


then  t/(s)a0  and  the  crack  arrests.  Given  values  of  iQ  and  Q, 
numerical  solutions  of  (5.5)  and  (5.9)  can  be  obtained,  with  the  value  of 
sr  being  determined  from  (4.16).  The  results  of  these  investigations 
have  shown  that  crack  arrest  does  indeed  occur  in  the  time  interval 
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( sr »sd ) *  However  <^ue  to  the  complicated  nature  of  I(s),  it  is  unclear 
whether  arrest  always  takes  place  in  the  interval  (sr,sd),  particularly 
for  very  large  Q. 

If  in  the  original  formulation  of  the  problem,  the  variables  x,y,  s 

and  w  are  made  dimensionless  by  considering  instead  x/2,  ,  y/2,  ,  s 1 1 

OOO 

and  w/wq,  then  the  governing  differential  equation  (4.1)  and  the  initial 
and  boundary  conditions  (3.3),  (4.2),  (4.3)  are  independent  of  all  para¬ 
meters  appearing  in  the  problem.  In  the  failure  criterion,  from  (2.3), 

(2.5)  and  (5.1),  the  only  variable  parameter  appearing  is  Q,  defined  in 
(5.3).  In  Fig.  3,  the  normalized  crack  tip  velocity  2/(s)  is  plotted 
against  s/lQ  for  Q  =*  1 . 5  and  crack  arrest  occurs  before  the  second 
reflected  wave  reaches  the  tip.  Similar  results  were  obtained  for  other 
values  of  Q,  both  larger  and  smaller  that  Q  ® 1.5. 

In  Fig.  4  the  dynamic  stress  intensity  factor  ( s ) ,  normalized 
by  its  Initial  static  value  given  in  (3.8),  is  shown  plotted  against 

s /lQ,  OstSiSSj,  again  for  Q  *  1 . 5 .  ( s )  drops  discontinuously  at 

s»0  from  Its  Initial  static  value  and  then  increases  due  to  decreasing 
crack  velocity  up  to  arrest.  After  arrest  occurs,  the  stress  intensity 
factor  decreases  until  the  arrival  of  the  next  reflected  stress  wave 
at  s  *  s^  .  The  effect  of  this  stress  wave  arriving  at  s » sd  is 
to  reduce  the  rate  at  which  the  stress  intensity  factor  is  decreasing  and 
possibly  even  to  begin  increasing  it  after  some  time.  The  effects  of  the 
succeeding  reflected  stress  waves  upon  the  stress  intensity  factor  become 
increasingly  more  difficult  to  analyze  by  this  method.  We  have  not  attempted 
such  an  analysis. 


THE  NUMERICAL  SOLUTION  OF  THE  EDGE  CRACK  PROBLEM 


6.  Introduction 

In  this  section  we  discuss  some  of  the  techniques  which  have  been 
used  as  well  as  the  difficulties  associated  with  the  numerical  solution 
of  rapidly  moving  crack  problems.  For  given  initial  and  boundary  condi¬ 
tions,  numerical  studies  of  dynamic  crack  propagation  problems  are 
generally  conducted  in  one  of  two  different  ways.  First,  given  the  crack 
tip  position  or  velocity  history  which  may  have  been  experimentally  ob¬ 
served  or  is  simply  prescribed,  then  by  forcing  the  numerical  model  to 
respond  in  exactly  the  same  manner,  one  can  calculate  items  of  interest 
such  as  the  dynamic  stress  intensity  factor  at  the  tip  or  possibly  use 
the  model  to  study  material  properties  from  laboratory  tests.  Alternatively, 
given  a  failure  criterion  as  In  (5.1),  one  can  compute  the  resulting  crack 
growth  history  and  the  crack  arrest  point  if  arrest  occurs.  This  may  be 
used  to  study  the  integrity  of  a  structural  component. 

Most  of  the  existing  numerical  work  done  in  simulating  two-dimensional 
crack  propagation  In  linear  elastic  bodies  has  been  performed  by  either 
the  finite  difference  or  the  finite  element  method.  In  a  critical  review 
of  numerical  solution  techniques  in  dynamic  fracture  mechanics  [14]  in 
1978,  Kanninen  discussed  some  of  the  deficiencies  inherent  in  both  methods 
and  stated  that  "Currently,  finite  elements  are  better  suited  for  initiation 
of  growth  calculations  while  finite  differences  can  better  cope  with  crack 
propagation  and  arrest".  Since  then  further  contributions  have  been  made 
to  both  of  these  methods  with  particular  emphasis  on  representing  the 
motion  of  the  crack  tip  and  the  singular  stress  field  locally  near  the 
tip  as  accurately  as  possible.  We  now  briefly  review  some  of  the  work 


which  has  been  done. 

In  the  finite  element  method,  the  modeling  of  crack  propagation 
problems  using  conventional  elements  everywhere  has  been  found  to  be  not 
sufficiently  accurate  [14].  Taking  advantage  of  special  elements  near  the 
tip  which  incorporate  the  characteristic  linear  elastic  crack  tip  singular¬ 
ity  in  the  stress  field,  the  finite  element  method  has  been  applied  suc¬ 
cessfully  to  equilibrium  crack  problems  [15]  and  to  crack  growth  initiation 
problems  [16].  Incorporating  the  singular  element  idea  at  a  moving  tip  in 
a  propagating  crack  problem  has  been  somewhat  more  difficult.  The  particu¬ 
lar  choice  of  singular  element  used  can  be  very  important  and  also  the 
advancement  of  the  tip  by  sudden  node  release  can  be  a  source  of  large 
transient  oscillations  in  the  numerical  model  [17,18]. 

Aoki  et  al  [19,20]  use  a  singular  element  incorporating  the  lowest 
order  asymptotic  displacement  eigenfunction  at  the  moving  crack  tip  and 
as  the  tip  gets  close  to  the  boundary  of  this  singular  element,  the  entire 
element  is  shifted  as  a  rigid  body  to  a  new  location.  Nishioka  and  Atluri 
[21,22]  use  a  linear  combination  of  crack  tip  displacement  eigenfunctions 
valid  for  a  constant  velocity  tip  to  construct  a  singular  crack  tip  element 
and  In  the  numerical  scheme  this  element  moves  with  the  tip  at  each  time 
step  being  surrounded  by  deforming  regular  finite  elements.  The  details 
Involved  however  are  somewhat  elaborate.  The  above  two  sets  of  authors  apply 
their  numerical  schemes  to  problems  involving  prescribed  crack  velocity 
[19,21]  and  also  to  cases  where  a  failure  criterion  is  given  and  the  crack 
velocity  history  Is  to  be  determined  [20,22].  Their  results  correlate 
reasonably  well  with  analytical  solutions  where  such  a  comparison  can  be 


The  finite  difference-  scheme  models  crack  growth  by  expressing  the 
elastodynamic  equations  of  motion  in  difference  form  at  each  node  point 
of  a  finite  grid  covering  the  body.  In  the  earlier  models  [23,24,25], 
crack  growth  is  simulated  by  intermittent  crack  extensions  of  magnitude 
one  mesh  size  which  implies  the  solution  of  a  sequence  of  transient  prob¬ 
lems,  each  contributing  to  a  disturbing  numerical  noise  with  no  physical 
counterpart.  In  [25]  a  transition  zone  is  used  at  the  moving  tip  to  aid 
in  remedying  the  above  difficulty.  Because  of  the  discontinuous  nature 
of  the  crack  motion  and  the  inadequate  representation  of  the  local  stress 
and  displacement  fields  near  the  moving  tip,  the  degree  of  resolution  and 
accuracy  obtained  in  the  numerical  results  from  these  models  is  not  suf¬ 
ficiently  good  for  extensive  investigations. 

To  improve  the  local  accuracy  at  the  moving  tip  and  in  particular  to 
overcome  the  problem  of  intermittent  crack  extension,  Aboudi[26]  trans¬ 
forms  the  equations  of  motion  to  a  moving  coordinate  system  with  origin 
at  the  moving  crack  tip.  A  finite  difference  scheme  involving  the  crack 
tip  velocity  and  acceleration  is  developed  for  both  the  plane  strain  and 
anti-plane  shear  cases  and  is  applied  to  a  number  of  problems  with  pre¬ 
scribed  crack  history.  Shmuely  and  Perl  [27,28]  introduce  in  the  vicinity 
of  the  propagating  crack  tip  a  local  moving  grid  which  moves  with  the 
velocity  of  the  tip  and  has  its  origin  of  coordinates  at  the  moving  tip. 
Outside  of  the  moving  grid,  a  stationary  grid  is  used  which  readily 
accommodates  the  boundary  conditions  for  finite  specimens.  Separate  dif¬ 
ference  schemes  are  established  governing  the  motion  in  both  the  moving 
and  stationary  domains.  As  In  the  model  by  Aboudi ,  the  difference  scheme 
for  the  moving  domain  directly  involves  the  crack  tip  velocity  and 


acceleration,  and  the  method  Is  applied  to  a  number  of  problems  also  with 
prescribed  crack  history.  In  the  prediction  case  of  crack  growth  embodying 
a  failure  criterion,  the  objective  is  the  determination  of  the  crack  tip 
position  and  velocity.  Finding  the  crack  tip  acceleration  is  an  order  of 
magnitude  more  difficult  and  requires  that  the  velocity  of  the  tip  be 
evaluated  with  a  high  degree  of  accuracy  and  resolution.  Hence  because 
the  acceleration  and  velocity  of  the  tip  occur  explicitly  in  both  of  the 
above  difference  schemes,  they  are  not  totally  suitable  for  the  determina¬ 
tion  of  crack  growth  history. 

Unlike  the  situation  in  the  finite  element  case,  little  effort  has 
been  made  in  finite  difference  simulations  to  incorporate  the  asymptotic 
nature  of  the  stress  and  displacement  fields  near  the  moving  crack  tip. 

In  not  doing  so,  the  singular  stress  field  near  the  tip  is  smoothed  out, 
thereby  affecting  the  numerical  solution  and  special  care  is  needed  in 
finding  the  stress  intensity  factor.  Burgers  [29]  uses  one  term  asymptotic 
expansions  about  the  crack  tip  position  to  approximate  the  stresses  and 
particle  velocity  in  the  region  near  the  tip  of  a  propagating  crack  in 
anti-plane  shear.  The  difference  equations  are  given  in  terms  of  the 
three  displacement  gradients  and  are  obtained  on  integrating  the  equations 
of  motion  along  bicharacteristic  strips.  The  stress  intensity  factor  is 
computed  for  both  steady  and  unsteady  prescribed  crack  tip  velocities  and 
compared  with  known  values.  The  results  are  reasonably  accurate  for  steady 
crack  propagation  but  fairly  large  errors  occur  for  the  unsteady  propaga¬ 
tion  cases. 

In  the  following  sections,  we  model  dynamic  crack  propagation  under 
anti-plane  shear  using  an  explicit  finite  difference  scheme  derived  from 


the  displacement  equation  of  motion.  To  improve  the  local  accuracy  near 
the  moving  crack  tip  at  which  singular  stresses  occur  and  to  smoothly 
extend  the  crack  at  each  time  step,  the  lowest  order  asymptotic  displace¬ 
ment  term  near  the  tip  is  used  in  conjunction  with  the  finite  difference 
scheme  to  simulate  the  problem  as  accurately  as  possible.  The  procedure 
is  applied  to  the  edge  crack  problem  discussed  in  Sections  2,3  and  4, 
finding  the  crack  tip  velocity  and  the  dynamic  stress  intensity  factor 
numerically.  To  check  the  accuracy  of  the  method, the  numerical  and 
analytical  results  are  compared. 


34 


This  section  deals  with  determining  the  behavior  of  the  out  of 
plane  displacement  w,  satisfying  (4.1)  and  (4.2),  near  the  tip  of  *n 
extending  crack  in  anti-plane  shear  and  employs  a  technique  first  used 
In  elastostatlcs  by  Kneln  [30]  and  Williams  [31].  A  similar  approach 
was  utilized  by  Freund  and  Clifton  [32]  in  analyzing  the  near-tip  dis¬ 
placement  field  for  non-uniform  crack  propagation  in  plane  strain.  The 
procedure  used  here  follows  that  of  Achenbach  and  Bazant  [33]  who  consider 
both  the  plane  and  anti-plane  situations  of  non-uniform  crack  extension 
and  obtain  the  lowest  order  near-tip  displacement  fields.  The  second 
order  correction  term  in  the  displacement  field  near  the  crack  tip  is 
established  here  for  non-uniform  crack  motion  in  anti-plane  shear.  In 
Implementing  the  numerical  simulation,  we  will  use  this  more  refined 
approximation  during  the  initial  stage  of  motion  of  the  crack. 

Consider  the  anti-plane  problem  of  a  straight  crack  extending  non- 
uniformly  in  Its  plane  along  the  x-axIs  say,  with  the  crack  tip  position 
given  by  x*Jl(s)  at  time  t,  where  as  before  we  let  s=ct,  c  being 
the  shear  wave  speed.  Introducing  the  local  polar  coordinate  system 
(R,0),  -  its  6  sir  at  the  moving  crack  tip  (Fig.  5)  through 

x-t(s)  ■  R  cos  9  ,  y*Rsin9  ,  (7.1) 

the  displacement  field  w*w(R,0,s)  from  (4.1)  satisfies  the  differential 
equation; 


wv»v"»vtv «’v»v« ■;  «y> *_•  v-  v *  • 


r-jrw  -  -w-g-w-  *  t  —  3  — Ti  a  -  •  -  v  ’  J 
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p¥^>R  f<^i)w,e  ■  2<*'  (,)«.,  e)  w,Rs  ♦  2(tLUh!m)w, 
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-  U'(s)cos  9  )w,R  +  ^-(s-^sin9)  W,( 


It  is  now  assumed  that 


w  %  Rmf(9,s)  as  R  -♦■  0,  -irs0iTT,  s^O 


(7.3) 


for  some  non-negative  exponent  m.  We  seek  the  smallest  value  of  m^O 
for  which  the  above  assumption  is  asymptotically  consistent  with  the 
differential  equation  (7.2)  and  the  traction-free  boundary  condition  (4.2) 
on  the  crack  surface.  Calculating  the  derivatives  of  w(R,9,s)  given  in 
(7.3)  asymptotically  as  R-*-0,  inserting  the  results  into  (7.2)  and  setting 
the  dominant  term  in  R  equal  to  zero,  one  finds  the  following  ordinary 
differential  equation  for  f(9,s): 


P(e)f,00-  (m-l)P,0(e)f,0  +  {m(m-l)^^Vm(2-m)P(e)}f -0  ,  (7.4) 

where  ^(9)  and  8  are  given  by 


( 


2  v 

l^)cos20,  8s6(s)s(l-i'(s)^  . 


P(9)  -  1  + 


(7.5) 
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To  solve  (7.4),  a  new  variable  0  is  introduced  through 

tan  0=8  tan  0  .  (7.6) 

For  8  given  in  (7.5)  and  0^Jl'(s)<l,  0  varies  continuously  from 
-ir  to  1  as  0  varies  from  -n  to  n.  In  terms  of  ©,  "F(e)  now 
becomes 

P(0)  »F(8(0))  =  (1  -  (1-B2)cos2©)-1  .  (7.7) 

Making  the  substitution 

m 

f(0,s)  *  P(0)7V(0,s)  (7.8) 

in  (7.4)  gives  the  following  differential  equation  for  v; 

V>6e  +  m2V  -  0  .  (7.9) 

The  traction-free  boundary  condition  on  the  crack  surface  gives 

f,.(±TT,s)*0  which  implies  V ,_ ( ±tt , s)  *  0  .  (7.10) 

For  displacement  fields  anti-symmetric  about  y  =  0,  as  will  be  considered 
here,  one  has  w*0  on  y»0  giving  V(0,s)  =  0.  The  appropriate  V 
satisfying  (7.9)  is  then 

V(©,s)  *  A(s)sinm©  ,  (7.11) 

A(s)  being  an  arbitrary  function  of  s  and  t(s).  The  smallest  value 
of  m  in  (7.11)  satisfying  (7.10)  is  m*J  .  Combining  (7. 3), (7. 8) 
and  (7.11)  now  gives  (see  also  [33]), 


w  ^  A(s)R^P(©)^sin  i©  as  R-*-0.  (7.12) 

To  find  the  second  order  correction  term,  let 

w  n.  R*f(0,s)  +  Rng(0,s)  as  R  +  0,  n  > \  ,  (7.13) 

where  f(0,s)  is  given  by  (7.8),  (7.11)  and  g(0,s)  and  n  are  to  be 
determined.  Again  obtaining  the  derivatives  of  w  asymptotically  as 
R->-0  from  (7.13),  putting  them  into  (7.2)  and  considering  the  dominant 
term  only,  it  can  be  shown  that  no  non-trival  function  g(0,s)  exists  for 
l/2<n<3/2.  For  n  =  3/2,  setting  the  dominant  term  in  R  to  zero  gives 

?( 9 ) g.00- ]■ F* 9 ( 6 } 9 * e +  9 >  39  a 

-  \  {2A'(sH'(s)  +  A(s)£'(s)){^(cos  0)f  -  (sin  0)f,0}  (7.14) 

8 

-  A(s  )l*  (s){j(cos  0)f-(sin0)f,0]»s 
8 

Introducing  @  and  P(©)  as  given  by  (7.6)  and  (7.7)  respectively  and 
making  the  substitution 

g(0,s)  *  P(©)5  Z(©,s)  ,  (7.15) 


equation  (7.14)  reduces  to 


Z,©0  + 


jZ-  ^  [2A'(s)i'(s)t/^jA(s)t'(s)]s(nj 


i^O.A(s)t'(s)s1n  |  6  .  (7.16) 
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The  boundary  conditions  are  Z,  (±tt,s)  =  0  and  Z(0,s)  =0  as  before. 

© 

Solving  (7.16)  now  for  Z(0,s)  and  using  (7.13),  (7.15)  gives  the  following 
expression  for  w  in  terms  of  0  and  s  as  R-*-Q; 

w  ^  A(s)R^P(0)^sin  0 

+  R^P(e)*^~  A(s)i*(s)^ sin  j-s  +  B(s)siny© 

+  illSj.  A(sU'(s)s1nl  el  .  (7.17) 
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Here  A(s),  B(s)  are  as  yet  arbitrary  functions  of  s,  and  t(s)  also 
remains  unspecified.  The  determination  of  A(s)  and  B(s)  depends  upon 
the  region  involved  and  the  particular  loading  on  the  body.  The  error  in 
(7.17)  is  0(R5/2)  as  R-0. 

Using  the  first  term  In  (7.17)  and  the  definition  of  the  dynamic 
anti-plane  shear  stress  Intensity  factor  ( s )  given  in  (2.6),  it 

can  be  established  that 

Kd(s)  -  u(£)Ma(s)  .  (7.18) 

Then  the  expression  for  the  energy  release  rate  given  in  (2.7)  becomes 

G(s)»^U(s)2  .  (7.19) 

By  virtue  of  the  failure  criterion  (5.1)  and  (7.19),  A(s)  is  determined 
explicitly  and  is  given  by 
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A(s)  “  for  a11  s’  (7.20) 

where  r  is  the  specific  fracture  energy,  assumed  constant.  A ( s )  is 
then  constant  for  all  s,  being  independent  of  2,(s)  also. 

For  the  anti-plane  shear  edge  crack  problem  discussed  in  Sections  3 
and  4,  A(s)  becomes 


o  rn~  W 

A  =  A(s)  =  ilf - t-  for  all  s,  (7.21) 

ttQ2  £2 

where  Q  is  the  dimensionless  parameter  given  by  (5.3),  %  is  the 

initial  crack  length  and  wQ  is  the  fixed  end  displacement  at  x  =  0  for 
y>0  (Fig.  1(a)). 

In  a  purely  static  problem,  for  a  straight  crack  along  the  x  axis 
with  local  polar  coordinates  r  and  a,  -ttscistt  situated  at  the  fixed 
crack  tip,  the  local  displacement  field  at  the  crack  tip  under  conditions 
of  anti-plane  shear  is 


w  ■  D^r^sin^-ct  +  Dgr^sin  j  a  +  0(r^72)  as  r-*-0 


(7.22) 


This  result  may  be  established  easily  by  employing  the  approach  used  in 
the  moving  crack  tip  case.  Considering  the  initial  static  displacement 
field  (3.3)  for  the  edge  crack  problem  in  Section  3,  it  can  easily  be 
shown  that 

~  "  w_ 
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(7.23) 


40 

Across  the  cylindrical  wavefront  r  =  s  generated  by  the  initial 
motion  of  the  crack  tip,  where  r  is  the  radial  distance  from  the  initial 
position  of  the  crack  tip  at  x  =  £o,  the  following  two  conditions  must 
hold; 


(i)  w  is  continuous  across  r  =  s  (7.2 

(ii)r^  +  j-  is  continuous  across  r = s,  (7.2 

where  ji  is  the  unit  outward  normal  to  r  =  s.  Since  the  material  must 
maintain  its  physical  integrity  at  the  wavefront,  equation  (7.24)  is 
necessary.  Equation  (7.25)  guarantees  that  the  impulse-momentum  relation 
at  the  wavefront  is  satisfied.  If  (7.24)  holds  across  r*s,  it  can 
easily  be  verified  that  (7.25)  is  also  valid. 

Equation  (7.24)  for  s«  1  can  be  used  to  find  fc'(0),  2*(0)  and 
B(0).  Equating  (7.17),  (7.21)  with  (7.22),  (7.23)  on  r*s  correct  to 
O(s^)  for  s«  1  gives 

i'( 0)  -  .  (7.2 

(r+i 

This  is  in  agreement  with  equation  (5.6)  which  results  from  the  exact 

differential  equation  governing  the  motion  of  the  tip.  Matching  correct 
3/2 

to  0(s  )  on  r*s  for  s«  1  gives 

Z'(0)  *  -  l*  (0)(1  -  £/(0)2)/4£q  (7.2 

and 

B  ( 0 )  »  A[16  +  5Jl'(0)  +  3H'(0)(1  -  l*  (0)2)]/l92  lQ  (1  -  ^(O)2)  , 


(7.2 


8.  Finite  difference  formulation  of  the  edge  crack  problem 


Here  we  consider  the  numerical  solution  of  the  dynamic  anti-plane 
shear  edge  crack  problem  formulated  and  discussed  in  Sections  2,3  and  4. 

The  two-dimensional  wave  equation  (4.1)  governing  the  motion  of  a  body  in 
a  state  of  infinitesimal  anti-plane  shear  must  be  satisfied  in  the  half¬ 
plane  xiO,  excluding  the  crack.  The  body  is  assumed  to  be  initially  at 
rest  and  the  boundary  conditions  are  given  by  (4.2)  and  (4.3).  The  initial 
static  solution  for  the  edge  crack  problem  is  given  in  (3.3).  The  problem 
as  formulated  is  anti-synmietric  about  y  =  0  and  hence  we  need  only  con¬ 
sider  the  region  x^O,  yiO. 

A.  Finite  difference  equations 


To  approximate  the  equation  of  motion  (4.1)  by  a  finite  difference 
equation  In  the  region  x*0,  yaO,  we  Introduce  a  rectangular  grid  of 
mesh  size  Ax*Ayxh  in  the  x  and  y*  directions  together  with  a  time 
Increment  As.  The  grid  nodes  are  located  at  (see  Fig.  6) 

x»(l-l)h,  yx(j-y)h;  i ,  j  =1,2,3, -  .  (8.1 

The  displacement  w  at  the  time  step  s*  k  As  and  at  a  mesh  point  (x,y) 
given  in  (8.1)  Is  denoted  by 

l  1  * j  *  1 ,2,3,---- 

w?  ,  •  w(x,y,s)  ;  (8.2 

k  =  0,1,2, - 

Let  P  denote  the  triad  ((1-l)h,  (j-j)h,  k  As  ).  Using  central  dif¬ 
ferences,  we  have 


Note  that  In  wf,j  the  comma  preceding  the  subscript  j  does  not  denote 
partial  differentiation. 
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For  w,  (P)+w,  (P),  we  use  the  nine  node  approximation: 

xx  yy 


.X,) 


w*xx(P)+w’yy(P)  "  C(w1-1,j-l  +  Vl,j+1  +  wi+1,j-l  +  Wi+l,j+l 

+wi.j+i  +wi+i »j*  -  (8-4) 


Equation  (8.4)  was  found  to  be  more  satisfactory  than  that  which  makes 
use  of  the  usual  central  difference  approximations  for  w,xx(P)  + w*yy(p) • 
Accordingly  from  (8.3)  and  (8.4),  we  get  the  following  finite  difference 
approximation  for  the  wave  equation  (4.1); 


k+1  _  ,_.k  _  Mk-1  A  p_  rfJ£  +  uk  . +wk  *  .  . +wk.,  . .,) 


'1J 


2w 


fj- Crr  cK-i.m  +  i-i  ,j+i +  wi+i  »j-i  + Vi,j+i 


+  4(wi-1,j+wi,j-l 


1  I 


(8.5) 


where  p*As/h  .  Equation  (8.5) is  a  three-level  explicit  difference 
scheme,  accurate  to  0(h^  +as^),  and  hence  it  is  possible  to  compute  the  dis¬ 
placement  at  time  step  s  =  (k  +  l)As  whenever  its  values  at  the  two  previous 
time  steps  s  *  (k  -  1)as  and  s  *  k  as  are  known  everywhere. 

The  boundary  condition  (4.3),  that  is  w  *  wQ  on  x  *  0,  ysO,  gives 


j  »  1,2,3,-— 
k  -  0,1,2,-— 


(8.6) 
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To  accommodate  the  boundary  conditions  on  ys0,  a  line  of  auxiliary  grid 


points  located  at  y»-h/2  outside  the  medium  is  introduced  (Fig.  6). 
The  crack  plane  is  on  y*0  and  the  crack  length  at  time  s  is  i{s) 
as  before.  For  Oaex<l(s)  the  free  surface  boundary  condition  yields 


w,y(x,0,s)  *  0,  which  on  approximating  by  central  differences  across  y-0 
gives 


w 


(8.7) 


for  all  i  such  that  0 s  (i-1  )h  <  2,(k  as ) .  For  x2l(s),  anti -symmetry  with 
respect  to  y*0  gives  w(x,0,s)=0  and  also  w(x,-h/2,s)  * -w(x,h/2 ,s) , 


which  yields 

wi,-l  a  ~wi  ,1  *  (8*8) 

for  all  i  such  that  (i-1  )h  >  £( k  As)  . 

Placing  the  grid  points  on  y»-h/2,  h/2,  3h/2, - was  foujid  superio 

to  putting  them  on  y*0,  h,  2h, - ,  thereby  having  grid  nodes  on  the  crack 

plane  and  in  the  path  of  the  crack  tip.  Equations  (8.7)  and  (8.8)  approxi¬ 
mate  the  boundary  conditions  on  y*0  for  x<£(s)  and  xii(s)  re- 
spectively,  correct  to  0(h)  away  from  the  crack  tip.  Because  of  the 
singular  displacement  gradients  at  the  tip  Itself,  the  finite  difference 
equation  (8.5)  and  the  corresponding  boundary  conditions  (8.7)  and  (8.8) 
cannot  be  expected  to  be  very  accurate  close  to  the  tip.  Additional 
measures  are  necessary  to  accurately  approximate  the  displacement  field 
near  the  moving  crack  tip. 


V£| 


i 
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To  apply  the  proposed  finite  difference  formula  (8.5),  stability 
criteria  must  be  established  which  guarantee  that  the  difference  between 
the  numerical  and  theoretical  solution  of  the  difference  equation  (8.5) 
does  not  grow  with  time.  The  analysis  performed  here  follows  von  Neumann's 
method  [34],  As  will  be  shown  in  the  following,  the  stability  requirement 
will  impose  restrictions  on  the  possible  values  of  the  time  increment  As 
in  relation  to  the  mesh  size  h  and  possibly  other  variables  in  the 
problem. 

The  stability  of  (8.5)  may  be  examined  by  assuming  that  an  error 

z?  .  exists  at  each  mesh  point  (k*0,  1,2, - ).  Since  the  difference 

* 

equation  (8.5)  is  linear  and  has  constant  coefficients,  the  growthof  these 

Ic  k 

errors  Is  governed  by  (8.5)  with  w^  ^  replaced  by  j  •  Performing 
a  harmonic  decomposition  of  the  errors  at  a  given  time  level,  say  s*0. 


gives 


r®  a  V  A  a 

MJ  //ij6 


(-l)*B,,.ih  (-l)*y.  Jh 


(8.9) 


where  the  frequencies  B^j  and  y^.  are  arbitrary  in  general  and  the 
summation  over  1  and  j  ranges  over  all  the  grid  points  (finite  in 
number)  In  the  x  and  y  directions,  respectively.  Since  the  difference 
equation  is  linear,  to  examine  the  growth  of  these  errors,  it  is  necessary 


to  consider  only  the  single  term 


(-l^Bih  (-l)Jyjh 


(8. 1C 


where  B  and  y  can  be  any  real  numbers. 
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To  Investigate  the  error  propagation,  we  now  seek  a  solution  to 
k  k 

(8.5)  with  w.  .  replaced  by  z.  .  which  reduces  to  (8.10)  at  s  =  0. 
• » J  • »J 

A  solution  of  the  difference  system  (8.5),  obtained  by  separating  the 
variables  with  (8.10)  in  mind,  is 


z*  .  =  R,e 
i.J  k 


(-l)*Bih  (-l)iyjh 


(8.11) 


where  satisfies 


-  2BRk  +  R^_1  *  0  ,  B-l-pfcA/2  , 


(8.12) 


A  «  4[sin2(^)-|-s1n2<^)sin2(^)+sin2(^)]  .  (8.13) 


The  general  solution  to  (8.12)  is  found  on  letting  R^n  which  gives 


ry  -  2Brj  +  1  *  0  , 


(8.14) 


or  n  *  B  ±  VS^T  .  The  errors  z*  ^  will  not  grow  provided  jn|  <1  , 


which  leads  to 


0  <  pcA  <  4  . 


(8.15) 


From  (8.13),  It  may  be  easily  verified  that  0^  As  16/3  ,  which  on  using 
(8.15)  gives  p</J/2  .  Hence  round-off  error  will  not  grow  and  the  dif¬ 
ference  scheme  (8.5)  will  be  stable  provided 


As  /T 

p  T  T 


(8.16) 


/•  ,n  V  .  % 


9.  Numerical  simulation  procedure  of  crack  motion  and  results 


Simulation  procedure  and  crack  tip  velocity  determination 


At  each  time  step  during  the  numerical  process,  the  velocity  of 
the  crack  tip  must  be  determined.  The  crack  tip  velocity  is  needed  in 
calculating  the  position  of  the  crack  tip  at  the  next  time  step  and  also 
is  required  in  utilizing  the  asymptotic  displacement  field  in  the  vicinity 
of  the  moving  crack  tip.  A  procedure  for  determining  the  crack  tip  veloc¬ 
ity  is  now  presented. 

Consider  a  region  0g  in  the  (x,y)  -  plane  containing  the  crack  tip  in 
its  interior  and  hence  having  a  portion  of  the  crack  surface  as  part  of  its 
boundary  (Fig.  5).  The  region  Ds  is  assumed  to  move  with  the  velocity 
of  the  tip  parallel  to  the  x-axis.  The  energy  contained  within  D$  at 
time  s  is  given  by 


E0  -  j  (W  ♦  K)dxdy  , 


(9.1) 


where  W  and  K  represent  the  potential  and  kinetic  energies  respectively 
of  the  body  per  unit  volume.  For  any  closed  contour  Cs  surrounding  the 
crack  tip  and  moving  with  the  velocity  of  the  tip  parallel  to  the  x-axis, 
the  rate  at  which  energy  is  entering  through  C$  is  given  by  the  contour 


integral 


Ec  *  c|  Cu|J|f  +  t'(s)(m-K)rgdCs  ,  (9.2) 


where  ii*(n  ,n  )  is  the  unit  outward  normal  to  C  and  c  is  the  shear 

*  y  * 

wave  speed.  Letting  C$  be  the  boundary  of  0$  excluding  the  crack 
surface,  then  the  flux  of  energy  through  C  given  in  (9.2)  must  equal 


the  rate  of  change  of  kinetic  and  potential  energy  within  D$  plus  the 

rate  at  which  energy  is  being  dissipated  by  the  moving  crack  tip.  Note 

that  there  is  no  contribution  to  Er  in  (9.2)  from  a  traction- free 

Ls 

crack  surface  along  the  x-axis.  Hence  we  have  the  following  relation  (see 
[10]  also); 

G*'(s)  -I  [u  !£!£  ♦  *'(s)(W  +  K)nx]dCs  I  (w  +  K)dxdy  ,  (9.3) 

Cs  °s 


where  6  is  the  energy  release  rate.  Equation  (9.3)  may  be  rewritten  as 


|  (w  +  Kjn^C^i' 


(*>M(£l?)dCs-ai- 

Cs 


(W  +  K)dxdy  ,  (9.4) 


giving  an  equation  for  the  crack  tip  velocity  i'(s). 

For  a  linear  elastic  body  in  anti-plane  shear, 

W  *  7  (w*x  +  w.y)  and  K  -  c2  j  w,2  .  (9.5) 

The  energy  release  rate  G  in  (9.4)  is  given  by  (7.19)  and  (7.21).  Then 
knowing  the  gradients  of  the  displacement  in  a  region  about  the  moving 
crack  tip,  equation  (9.4)  enables  us  to  calculate  the  velocity  of  the  tip. 
This  must  be  done  using  the  numerical  displacement  field. 

As  was  noted  in  Section  8,  the  finite  difference  equation  (8.5) 
cannot  be  expected  to  be  very  accurate  at  the  mesh  points  closest  to  the 
moving  crack  tip.  To  alleviate  this  difficulty  the  lowest  order  asymptotic 
displacement  term  (7.12),  with  A(s)  given  by  (7.21),  is  used  at  the  two 
nodal  points  In  the  region  nearest  to  the  crack  tip  i.e.  at  (LN,1)  and 
(LN+1,1),  where  1  »  LN  and  j3l  denote  the  mesh  point  nearest  to  the 
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tip  and  to  the  left  of  it.  To  smooth  the  transition  from  the  asymptotic 

to  the  numerical  values  of  the  displacement  as  the  tip  moves,  a  linear 

combination  of  both  the  asymptotic  and  numerical  displacement  values  is 

used  at  the  four  nodal  points  (LN-1,1),  (LN,2),  (LN+1,2)  and  (LN+2,1) 

surrounding  the  tip.  The  linear  combination  used  is  [aw1?  .+(h-a)w  ]/h, 

1 ,  .1  as 

where  w..  denotes  the  lowest  order  near-tip  asymptotic  displacement 
given  in  (7.12)  and  a  s  [S,(  kAs) -(LN  -  1  )h]  is  the  distance  of  the  tip 
from  i  *  LN . 

Before  calculating  the  crack  tip  velocity  for  the  present  or  kth 
time  step,  the  lowest  order  asymptotic  near-tip  displacement  (7.12)  is 
used  as  outlined  in  the  preceding  paragraph  together  with  the  crack  tip 
velocity  of  the  previous  or  (k-l)st  time  step.  After  determining  the 
velocity  of  the  tip  for  the  kth  time  step,  equation  (7.12)  is  then 
used  again  with  this  new  crack  tip  velocity  to  provide  the  near-tip 
asymptotic  displacement  for  the  current  time  step.  Knowing  the  crack  tip 

t  h 

position  £.(  k  As )  and  velocity  ^(kAs)  at  the  k  time  step,  the 
position  of  the  tip  at  the  next  time  step  is  given  by 

£((k  +  l)As)  *  H(kAs)  +Jl/(kAs)AS  ,  (9.6) 

2 

correct  to  0(As  )  . 

To  initiate  the  motion  of  the  crack  tip,  the  near-tip  asymptotic 
displacement  including  the  second-order  correction  term  given  by  equation 
(7.17)  is  used  In  conjunction  with  (7.21),  (7.27)  and  (7.28)  at  the  ap¬ 
propriate  mesh  points  near  the  tip  for  the  first  four  time  steps.  The 
crack  tip  position  and  velocity  used  with  equation  (7.17)  for  these  first 
four  time  steps  is  given  by 
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2(0)  +  l'( 0)s  +  J  l'( 0)s2 


2'(s)  =  i'(0)  +  2'(o)s  , 


where  2/(0)  and  2#(0)  are  given  in  (7.26)  and  (7.27)  respectively. 

To  facilitate  further  the  initiation  of  the  numerical  process,  equations 
(9.7)  and  (9.8)  are  used  to  give  the  crack  tip  position  and  velocity  for 

the  first  ten  time  steps.  This  is  justified  for  sufficiently  small  time 

3  2 

increments  As  since  the  error  in  (9.7)  and  (9.8)  is  0(s  )  and  0(s  ) 
respectively  for  s«  1  .  Thereafter  the  crack  tip  velocity  is  found 
numerically. 

To  determine  the  crack  tip  velocity  from  equation  (9.4),  the  time 
rate  of  change  of  the  kinetic  and  potential  energy  in  the  region  D$, 
i.e.  the  last  term  in  (9.4),  must  be  evaluated  numerically.  However, 
because  of  transient  disturbances  introduced  into  the  numerical  values 
of  the  displacement  by  the  motion  of  the  crack  tip,  the  accurate  calcula¬ 
tion  of  the  last  term  in  (9.4)  proved  to  be  very  difficult.  Using  the 
lowest  order  asymptotic  expression  for  the  displacement  near  the  crack  tip 
(7.12),  it  can  be  shown  that 


(W  +  K)dxdy  *  0(h) 


as  h  ->■  0 


Then  taking  a  very  small  contour  Cs  about  the  crack  tip  and  neglecting 
the  last  term  in  (9.4)  gives  the  following  equation  which  approximates 
2'(s)  correct  to  0(h)  as  h->-0  , 


1  *  ■  V  ■ 


:,l».';r'f.».r,» — 3-; 
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(•-f  <W*,C)„xdC^.>  •*{<£&)  dC,  . 


(9.10) 


From  equation  (9.6)  the  crack  tip  position  for  the  next  time  can  then  be 
estimated  correct  to  0(h  ).  The  energy  release  rate  G  is  given  by 
(7.19),  (7.21)  and  substituting  for  W  and  K  from  (9.5),  equation 
(9.10)  becomes 


/2wo2  .  I  x2  ,9w  ,2  .  ,9w  v21  ,r 

yTTQ2-o  ^  9x  9y  ^  9s  nxd  y  ^ 


(sM(fS§)dcs 


(9.11) 


with  Q  given  In  (5.3). 

In  Fig.  6,  the  upper-half  of  the  contour  C$  used  in  finding  the 

velocity  of  the  crack  tip  is  illustrated.  The  contour  is  fixed  relative 

to  the  tip  for  all  time  since  it  moves  with  the  velocity  of  the  tip. 

From  the  symmetry  Inherent  in  the  problem,  it  is  sufficient  to  evaluate 

the  integrals  in  the  upper  half-plane.  The  vertical  segments  of  the 

contour  are  a  distance  2h  ahead  and  behind  the  crack  tip,  while  the 

horizontal  portion  is  along  ja3.  To  evaluate  the  contour  integrals  in 

(9.11),  the  displacement  gradients  w,  ,  w,  and  w,  must  be  found  at 

x  y  s 

points  along  the  contour.  At  the  mesh  points  adjacent  to  or  on  the 

contour,  the  spatial  gradients  w,  and  w,  are  found  using  the  central 

x  y 

difference  approximation.  To  find  w,g  at  the  grid  node  (i,j)  and  at 
the  time  step  s»kAs,  the  following  difference  approximation  is  used 


/  9w  \  k 

V3sM,j 


(9.12) 


The  error  in  (9.12)  is  O(As^)  and  for  the  (k-2)nc*  time  step  it  is 
only  necessary  to  keep  the  values  of  the  displacement  at  the  required 
mesh  points  near  the  tip.  The  integrals  along  the  top  segment  of  the 
contour  are  evaluated  using  the  three  point  Simpson's  rule,  while  the 
integrals  along  the  vertical  segments  are  evaluated  using  a  three-point 
formula  with  integration  points  at  y  =  h/2,  3h/2  and  5h/2,  which  is  exact 
for  cubic  polynomials.  At  a  requisite  point  on  the  contour,  the  displace¬ 
ment  gradients  are  found  using  a  linear  combination  of  the  values  of  the 
corresponding  gradients  at  the  two  nearby  grid  nodes,  one  to  the  right 
and  one  to  the  left  of  the  point  under  consideration.  This  gives  the  dis¬ 
placement  gradients  at  the  integration  points  on  the  contour  correct  to 
0(h  ).  The  normalized  crack  tip  velocity  j/(s)  can  now  be  found  from 
equation  (9.11), 

B.  Stress  intensity  factor  determination 

To  find  the  dynamic  stress  intensity  factor  at  the  moving 
crack  tip,  we  note  that  from  equations  (7.18)and  (7.5)  that 

Kd(s)  »  (J)*u(l  (s)Z)*A(s)  ,  (9.13) 

where  A(s)  is  the  time-dependent  function  in  the  expression  for  the 
lowest  order  asymptotic  displacement  near  the  tip  in  (7.12).  Normalizing 
Kd(s)  with  respect  to  the  initial  static  stress  intensity  factor  K$t 
given  in  (3.8),  equation  (9.13)  becomes 


As  the  crack  is  propagating,  the  failure  criterion  (5.1)  applies;  it 
yields  the  previously  derived  expression  in  equation  (7.21)  for  A(s). 
Using  (7.21),  equation  (9.14)  reduces  to 


Kj(s)  -  ,  . 

4 - (1  -  ^  (s)Z)«  /  , 


(9.15 


and  this  equation  is  valid  while  the  crack  is  extending  i.e.  ?/(s)>0. 
Hence  up  to  the  time  of  arrest,  the  stress  intensity  factor  may  be  easily 
calculated,  once  the  crack  tip  velocity  is  determined  from  the  numerical 
solution. 

After  arrest  the  failure  criterion  is  no  longer  valid,  and  the  stress 

intensity  factor  must  be  evaluated  directly  from  the  numerical  solution. 

Equation  (7.12)  is  still  valid  after  arrest  with  i'(s)  equal  to  zero. 

Consider  now  the  upper-half  of  the  contour  C  ,  used  in  evaluating  the 
#  ^ 

crack  tip  velocity  l*  (s)  in  (9.11),  which  we  denote  by  C*  .  The  posi¬ 
tion  of  C*  is  now  fixed  after  arrest  relative  to  the  stationary  tip. 

To  evaluate  the  stress  intensity  factor,  we  consider  the  contour  integral 


Using  for  w  the  lowest  order  asymptotic  displacement  term  near  the  tip 
given  in  (7.12)  and  evaluating  the  contour  integral  analytically  gives 


.  wdC*  -  [1  -  (5+(41)M-5(5  +  (41)*)'*]A(s)h^  +  0(h*)  .  (! 


Rewriting  equation  (9.17)  gives  the  following  approximation  to  A(s); 

A(s)  =  —  h“*  [1  -  (5+{41)M  -  5(5+(41  )*)“*]  1  wdC*-  (9*18) 
4SZ  C* 

The  contour  integral  in  (9.18)  is  calculated  in  a  similar  manner  to  the 
way  in  which  the  contour  integrals  in  (9.11)  were  evaluated  when  determin¬ 
ing  the  crack  tip  velocity.  This  enables  us  to  determine  A(s).  Equation 
(9.14)  with  2,'(s)  set  equal  to  zero  and  equation  (9.18)  then  give  the 
stress  Intensity  factor  after  arrest. 

C.  Numerical  results 

In  the  numerical  results  presented  here,  the  exact  initial  static 
displacement  field  (3.3)  Is  used  at  all  grid  points  in  the  mesh  at  s*0. 

The  time  Increment  taken  at  each  time  step  is  As*h/2,  which  satisfies 
the  stability  requirement  (8.16).  Considering  the  dimensionless  variables 
x/Jl0,  y/*0,  s/lQ  and  w/wQ  as  done  in  Section  5  for  the  analytical  case,  then 
the  problem  depends  only  on  the  dimensionless  parameter  Q,  defined  in 
(5.3).  As  in  the  analytical  case,  Q  s  1 .5  is  taken  here. 

Plots  of  the  crack  tip  velocity  and  crack  length  evaluated  numerically 
up  to  the  time  of  arrest  are  shown  In  Figs.  7  and  8.  The  mesh  size  employed 
Is  Ax*Ay*h*.02,  while  the  number  of  grid  points  used  for  the  quarter- 
plane  xiO,  yiO  is  (150,90).  For  the  value  of  h  under  consideration 
this  ensures  that  the  outer  boundary  effects  arising  from  the  finite  grid 
do  not  affect  the  results.  In  Fig.  7  the  normalized  crack  tip  velocity 
Jt# (s)»i(t)/c  ,  calculated  numerically  using  (9.11),  is  plotted  against 
s * ct  and  compared  with  the  exact  tip  velocity  acquired  in  Section  5. 
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The  numerically  evaluated  velocity  tends  to  oscillate  si ightly  as  the 
crack  tip  traverses  the  grid,  but  overall  the  agreement  with  the  exact 
solution  is  quite  good.  The  reflected  wavefront  from  the  boundary  at 
x *  0  is  smoothed  out  in  the  numerical  solution  and  hence  the  tip  velocity 
as  determined  from  the  numerical  solution  does  not  exhibit  the  abrupt 
corner  present  in  the  exact  solution  when  this  wave  reaches  the  crack  tip. 
This  is  to  be  expected,  since  the  difference  sheme  cannot  propagate  a 
sharp  wavefront.  In  Fig.  8  the  normalized  crack  length  l(s)/iQ  is  plotted 
against  s/lQ~  ct/iQ  up  to  arrest  and  shows  excellent  agreement  with 
the  exact  solution. 

The  numerically  calculated  stress  intensity  factor  Kd(s),  normalized 
by  the  initial  static  stress  intensity  factor  Kst,  is  shown  in  Fig.  9  up 
to  the  time  that  the  second  reflected  wave  reaches  the  tip  from  the  end 
x*0.  The  mesh  size  used  in  Fig.  9  is  h*  .04  while  the  number  of  grid 
points  is  (150,90)  as  before.  Before  crack  arrest  Kd(s)  is  determined 
using  equation  (9.15)  together  with  the  numerically  evaluated  crack  tip 
velocity,  while  after  arrest  equations  (9.14)  and  (9.18)  are  used.  The  reflected 
stress  waves  return  to  the  tip  slightly  faster  in  the  numerical  than  in 
the  analytical  solution  because  the  wavefronts  are  being  dispersed  somewhat 
in  the  numerical  scheme  as  time  increases. 

Overall  the  numerical  results  agree  remarkably  well  with  the  exact 
solution.  The  accuracy  of  these  results  increases  for  smaller  values  of 
h.  Also  the  smaller  the  value  of  As  relative  to  h  ,  the  better  the 
results,  primarily  because  more  time  steps  are  then  used  in  traversing 
the  grid,  and  this  leads  to  a  smoother  velocity  profile.  Care  must  be 
exercised  in  using  the  asymptotic  solution  near  the  moving  tip  to  avoid 
introducing  unwanted  oscillations  into  the  crack  tip  velocity.  The  method 
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presented  herein  for  avoiding  such  oscillations  as  much  as  possible  was 
found  to  be  the  most  satisfactory. 
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10.  Further  appl ications 


All  of  the  work  presented  in  this  study  has  been  for  infinitesimal 
anti-plane  shear  deformations.  In  plane  strain  crack  propagation  problems, 
two  non-zero  displacement  components  are  present,  as  well  as  both  dilata¬ 
tion  and  shear  waves.  This,  together  with  more  complicated  wave  reflection 
patterns  at  the  specimen  boundaries  makes  plane  strain  problems,  especially 
for  finite  regions,  far  more  difficult.  Therefore  one  must  resort  to  a 
numerical  approach.  As  in  anti-plane  shear,  the  stress  and  particle 
velocity  fields  are  singular  at  the  moving  tip  and  hence  must  be  adequately 
approximated  in  a  numerical  simulation. 

The  lowest  order  asymptotic  displacement  field  near  the  tip  of  a 
propagating  crack  in  plane  strain  is  available  for  non-uniform  crack 
velocity  In  [32]  and  so  can  be  incorporated  into  a  suitable  finite  differ¬ 
ence  scheme  In  order  to  accurately  represent  the  elastodynamic  field  near 
a  propagating  crack  tip.  It  should  be  possible  to  accomplish  this  without 
a  lot  of  difficulty,  not  only  for  the  case  of  known  crack  history,  but 
where  the  crack  history  Is  to  be  determined  using  a  suitable  failure 
criterion  as  well. 

It  would  be  of  Interest  to  apply  the  finite  difference  technique 
Incorporating  the  local  displacement  field  near  a  crack  tip  to  some  plane 
strain  crack  propagation  problems  already  studied  and  compare  with  such 
exact  results  as  are  available  or  with  other  numerical  solutions.  An 
Interesting  problem  suitable  for  numerical  study  is  the  plane  strain 
analogue  of  the  edge  crack  problem  analyzed  herein.  For  a  constant  dis¬ 
placement  In  the  positive  y  direction  and  an  equal  and  opposite  displace¬ 
ment  in  the  negative  y  direction  at  the  boundary  x*0,  thus  separating 


the  crack  surfaces;  then  for  fixed  or  traction-free  conditions  normal  to 
the  end  xa0,  one  could  study  the  effect  of  the  individual  reflected 
stress  waves  from  x*0  on  the  motion  of  the  tip  and  check  if  and  when 
crack  arrest  might  occur. 

The  anti-plane  shear  work  discussed  in  this  investigation  has  been 
carried  out  for  infinitesimal  deformations  in  linear  elastic  materials. 
Knowles  [35]  has  studied  the  elastic  field  near  the  tip  of  a  stationary 
crack  for  a  class  of  fully  nonlinear  elastic  materials  subjected  to 
finite  deformations.  It  would  be  interesting  to  consider  the  dynamic 
problem  of  a  rapidly  extending  crack  in  a  nonlinear  elastic  material 
taking  finite  deformations  into  account.  The  local  displacement  field 
near  a  moving  crack  tip  may  be  derived,  as  in  Section  7,  for  different 
classes  of  nonlinear  materials.  This  could  then  be  used  with  a  finite 
difference  scheme  to  model  the  problem,  and  thus  to  study  the  effect  of 
nonlinearity  on  the  crack  tip  motion. 


REFERENCES 


[1]  B.V.  Kostrov,  Unsteady  propagation  of  longitudinal  shear  cracks. 

Journal  of  Applied  Mathematics  and  Mechanics  (English  translation 
of  PMfl),  30  (1966)  1241. 

[2]  J.O.  Eshelby,  The  elastic  field  of  a  crack  extending  non-uniformly 
under  general  anti -plane  loading.  Journal  of  the  Mechanics  of 
Physics  and  Solids,  17  (1969)  177. 

[3]  L.B.  Freund,  Crack  propagation  in  an  elastic  solid  subjected  to 
general  loading-11.  Non-uniform  rate  of  extension.  Journal  of  the 
Mechanics  of  Physics  and  Solids,  20  (1972)  141. 

[4]  B.V.  Kostrov,  Crack  propagation  at  variable  velocity.  International 
Journal  of  Fracture,  11  (1975)  47. 

[5]  T.C.  Kennedy  and  J.O.  Achenbach,  Extension  of  a  cut-induced  crack 
under  antiplane  loading.  Journal  of  Elasticity,  3  (1973)  277. 

[6]  J.O.  Achenbach,  Extension  of  a  crack  by  a  shear  wave,  ZAMP,  21  (1970) 

887.  - 

[7]  J.D.  Achenbach,  Crack  propagation  generated  by  a  horizontally 
polarized  shear  wave.  Journal  of  the  Mechanics  of  Physics  and  Solids, 

18  (1970)  947. 

[8]  A.E.  Molchanov  and  L.V.  Nikitin,  Dynamics  of  longitudinal  shear  cracks 
after  the  loss  of  stability.  Mechanics  of  Solids  (Enqlish  translation 
of  MTT),  7  (1972)  55. 

[9]  P.H.  Melville,  Effect  of  reflected  stress  waves  on  the  stress  intensity 
factor  of  an  anti-plane  strain  edge  crack  during  crack  propagation  and 
after  crack  arrest,  International  Journal  of  Fracture,  17(1981)  3. 

[10]  C.  Atkinson  and  J.D.  Eshelby,  The  flow  of  energy  into  the  tip  of  a 
moving  crack,  International  Journal  of  Fracture  ,4  (1968)  3. 

[11]  J.R.  Rice,  Mathematical  analysis  in  the  mechanics  of  fracture. 

Fracture,  Vol.  2,  edited  by  H.  Llebowitz,  Academic  Press.  New  York, 
(19W  191. 

[12]  J.D.  Achenbach,  Dynamic  effects  in  brittle  fracture.  Mechanics  Today, 
edited  by  S.  Nemat-Nasser,  Pergamon  Press,  New  York,  (1974)  1. 

[13]  J.  Mlklowltz,  The  Theory  of  Elastic  Waves  and  Waveguides,  North-Hoi  land. 
New  York,  (1978)  363. 


[14]  M.F.  Kannlnen,  A  critical  appraisal  of  solution  techniques  in 
dynamic  fracture  mechanics,  Numerical  Methods  in  Fracture  Mechanics, 
edited  by  A.R.  Luxmoore  and  D.R.J.  Owen,  University  College  of 
Swansea,  Wales,  (1978)  612. 

[15]  R.H.  Gallagher,  A  review  of  finite  element  techniques  in  fracture 
mechanics.  Numerical  Methods  in  Fracture  Mechanics,  edited  by 

A.R.  Luxmoore  and  D.R.J.  Owen,  University  College  of  Swansea,  Wales, 
(1978)  1. 

[16]  J.A.  Aberson,  J.M.  Anderson  and  W.W.  King,  Dynamic  analysis  of 
cracked  structures  using  singularity  finite  elements,  Elastodynamic 
Crack  Problems,  Mechanics  of  Fracture  4,  edited  by  G.C.  Sih,  Noordoff 
International  Publishing,  Leyden,  (1977)  249. 

[17]  J.A.  Aberson,  J.M.  Anderson  and  W.W.  King,  Singularity-element  simu¬ 
lation  of  crack  propagation.  Fast  Fracture  and  Crack  Arrest,  ASTM 
STP  627,  edited  by  G.T.  Hahn  and  M.F.  Kanninen,  American  Society  for 
Testing  and  Materials,  (1977)  123. 

[18]  W.W.  King,  J.F.  Malluck,  J.A.  Aberson  and  J.M.  Anderson,  Application 
of  running-crack  eigenfunctions  to  finite-element  simulation  of 
crack  propagation.  Mechanics  Research  Communications,  3  (1976)  197. 

[19]  S.  Aoki ,  K.  Klshimoto,  H.  Kondo  and  M.  Sakata,  Elastodynamic  analysis 
of  crack  by  finite  element  method  using  singular  element.  International 
Journal  of  Fracture,  14  (1978)  59. 

[20]  K.  Klshimoto,  S.  Aokt  and  M.  Sakata,  Computer  simulation  of  fast 
crack  propagation  in  brittle  material.  International  Journal  of 
Fracture,  16  (1980)  3. 

[21]  T.  Nlshloka  and  S.N.  Atlurl,  Numerical  modeling  dynamic  crack  propaga¬ 
tion  in  finite  bodies,  by  moving  singular  elements  —  Part  1 :  Formulation, 
Part  2:  Results,  Journal  of  Applied  Mechanics,  47  (1980)  570. 


[22]  T.  Nlshloka  and  S.N.  Atlurl,  Numerical  analysis  of  dynamic  crack 
propagation:  generation  and  prediction  studies.  Engineering  Fracture 
Mechanics,  16  (1982)  303. 

[23]  M.  Shmuely  and  D.  Peretz,  Static  and  dynamic  analysis  of  the  DCB 
problem  In  fracture  mechanics,  International  Journal  of  Solids  and 
Structures ,  12  (1976)  67. 

[24]  T.  Kanazawa  and  S.  Machida,  Energy  considerations  in  dynamic  crack 
propagation  and  arrest.  Proceedings ,  Fourth  International  Conference 
on  Fracture,  edited  by  O.M.R.  Taplin,  Waterloo,  Canada,  (1977)  223. 

[25]  H.  Stock 1  and  F.  Auer,  Oynamic  behavior  of  a  tensile  crack:  finite 
difference  simulation  of  fracture  experiments.  International  Journal 
of  Fracture,  12  (1976)  345. 


61 


[26]  J.  Aboudi ,  Numerical  solution  of  dynamic  stresses  induced  by  moving 
cracks.  Computer  Methods  in  Applied  Mechanics  and  Enqineerinq, 

9  (1976)~3T5T; - a - * 

[27]  M.  Shmuely  and  M.  Perl,  The  SMF2D  code  for  proper  simulation  of 
crack  propagation.  Crack  Arrest  Methodology  and  Applications,  ASTM 
STP  711,  edited  by  G.T.  Hahn  and  M.F.  Kanninen,  American  Society  for 
Testing  and  Materials,  (1980)  54. 

[28]  M.  Perl,  M.  Shmuely  and  D.  Peretz,  An  improved  dynamic  crack  propaga 

tion  simulation  in  the  SEN  specimen  by  the  SMF20  code.  International 
Journal  of  Fracture,  19  (1982)  17.  ~ 

[29]  P.  Burgers,  Dynamic  linear  elastic  crack  propagation  in  anti-plane 

shear  by  finite  differences.  International  Journal  of  Fracture,  16 
(1980)  261.  - 

[30]  M.  Knein,  Zur  theorie  dec  druckversuchs,  Abhandl ungen  aus  dem 
aerodynamischen  Institut  an  der  T.H.  Aachen,  1  (1927)  62. 

[31]  M.L.  Williams,  Stress  singularities  resulting  from  various  boundary 

conditions  in  angular  corners  in  extension.  Journal  of  Applied 
Mechanics,  19  (1952)  526.  - 

[32]  L.B.  Freund  and  R.J.  Clifton,  On  the  uniqueness  of  plane  elasto- 

dynamic  solutions  for  running  cracks.  Journal  of  Elasticity,  4 
(1974)  293.  - 

[33]  J.D.  Achenbach  and  Z.P.  Bazant,  Elastodynamic  near-tip  stress  and 
displacement  fields  for  rapidly  propagating  cracks  in  orthotropic 
materials.  Journal  of  Applied  Mechanics,  42  (1975)  183. 

[34]  A.R.  Mitchell  and  D.V.  Griffiths,  The  Finite  Difference  Method  in 
Partial  Differential  Equations,  Wiley,  (1980)  38. 

[35]  J.K.  Knowles,  The  finite  anti-plane  shear  field  near  the  tip  of 

a  crack  for  a  class  of  incompressible  elastic  solids,  international 
Journal  of  Fracture,  13  (1977)  611. 


V< i'v vv vv v*v  ■ •  %V- -.y.\nY -\v. 


0 


1.0 


3. 


2.0 

S/^o 


0 


Figure  3.  Normalized  crack  tip  velocity  fc'(s)*£(t)/c  as  a  function 
of  s/*0  a  ct/i,Q,  from  the  analytical  solution. 


*_•  V  V  V 


Unclassified 


S*Cy«'^v  C.-»SSi*>CAriON  Q»  this  P4G*  _  I 


|  REPORT  DOCUMENTATION  PAGE 

U  SSsiwHiTv  C^ASSIF.CAliON 

Unclassified 

lO.  «eSTRjCTlvl  MARKir.Ga 

2*.  SECURITY  CLASSIFICATION  iuThOHirv 

3.  O'STRIBLiTION.  AVAiuAdiL:  T  v  OF  <EaCRT 

2a.  OSCLASSlFiCATiQN.  DOWNGRADING  SChEOULE 

4.  PERFORMING  ORGANIZATION  REPORT  NVjMBSRISI 

5.  monitoring  organisation  A 6 PORT  \„\i8e»iS) 

6*  SAMS  OF  PERFORMING  ORGANIZATION 

California  Institute  of  Tech. 

ao.  OFFICE  SYMBOL 
■  /Tapp  Leapt*  i 

7»  NAME  OF  MONITORING  ORGANIZATION 

Office  of  Naval  Research 

6c.  AOQRCSS  iCil).  auu  and  /.IP  Cwd*i 

Pasadena,  California  91125 

7b.  AOOR6SS  /City,  ztatt  u/ia  ///*  Couci 

Washington,  D.C. 

feNAMl  OF  FUNOING/SFONSOMING 
ORGANIZATION 

80.  OFFICE  SYMBOL 

Ilf  appLccPIfl 

9.  PROCUREMENT  INSTRUMENT  iOENTiFiCaTion  NUMBER 

N00014-75-C-0196 

Sc.  AOQMESS  iCili,  olau  and  ZIF  Co dc< 

10.  SOURCE  O*  PUNOING  NOS. 

PROGRAM  j  PROJECT  i  *AS<  |  ’.VQRK  UNIT 

element  no.  j  no.  ,  no.  i  no. 

i  .  i 

11.  titlS  /afiuiM  Mcvmy  c/MH/ifgiioAi  TFltf  [J^bp^dUUn  Ciiii  a 

space  under  conditions  of  anti -plane  shear:  £n 

^rest  of  an  eppe  cracK  in  an  eiasxio  naiT- 
ilytical  and  humerical  results.  i 

I  12.  MRSONAL  AuThORISI 

| Timothy  C.  O'Sullivan 

1  13*  TYPE  OF  REPORT  130.  TIME  COVeRSO 

|  Technical  from  to 

14  OATE  OF  REPORT  .1  r  Mo..  Day)  15  PAGE  COUNT  j 

September  1983  68  j 

11  sumumintahv  NOTATION 

-  -  —  —  -  _ 

IT  CCSATI  COOCS 

18.  SU8j8CT  TC  AMS  •Continue  on  reverie  •(  ettery  end  identify  by  biocn  number # 

KEYWORDS:  crack  propagation,  crack  arrest,  dynamic  fracture. 

FifLO  1  GROUP  1  sue.  GR 

1  - 

i  ; 

!•.  ABSTRACT  'Continue  V«I  reverie  if  nrceiiery  itfinfi/y  by  bloc*  number* 


The  notion  of  an  edge  crack  extending  non-uniformly  in  an  elastic  half-space  under  condi-  1 
cions  of  anti-plane  shear  is  analyzed.  An  expression  for  the  stress  intensity  factor  I 

at  the  crack  tip  is  obtained,  and  an  energy  balance  crack  propagation  criterion  is  used  j 

to  find  the  equation  of  motion  of  Che  tip.  On  solving  this  equation  numerically,  it  is  I 

found  that  crack  arrest  occurs  before  the  second  reflected  wave  from  the  boundary  reaches  j 
the  tip.  | 

i 

In  the  second  half  of  this  investigation,  a  numerical  procedure  for  studying  anti-plane  ] 
shear  crack  propagation  problems  using  finite  differences  is  developed.  To  approximate  1 
Che  elastodynamic  field  as  accurately  as  possible  near  the  moving  crack  tip,  where  singu¬ 
lar  stresses  occur,  the  local  asymptotic  displacement  field  near  the  tip  is  incorporated 
into  the  finite  difference  scheme.  The  numerical  procedure  is  applied  to  the  edge  crack 
problem  analyzed  in  the  first  part  of  this  study,  and  the  numerical  and  exact  results  are 
compared. 

•  caTPi  ,;:  ;s  r  . .  j.r  ■  Ca  r  i ^ j 

3 


.9  «•  * 


•  •jj  r 


MS  -  ' 


^^■13  •  *« W •  V»  r i 


nor 


